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Abstract 

The Casimir free energy for a system of two dielectric concentric 
nonmagnetic spherical bodies is calculated with use of a quantum 
statistical mechanical method, at arbitrary temperature. By means 
of this rather novel method, which turns out to be quite powerful (we 
have shown this to be true in other situations also) , we consider first an 
explicit evaluation of the free energy for the static case, corresponding 
to zero Matsubara frequency (n = 0). Thereafter, the time-dependent 
case is examined. For comparison we consider the calculation of the 
free energy with use of the more commonly known field theoretical 
method, assuming for simplicity metallic boundary surfaces. 
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1 Introduction 



The Casimir problem for dielectrics - for general introductions see, for in- 
stance Refs. [1-3] - has turned out to be difficult to solve, in the presence of 
curved surfaces. The most typical example of a system of this sort is prob- 
ably that of a single nonmagnetic compact spherical ball, surrounded by a 
vacuum. ( Equivalent ly, one may imagine a spherical cavity in an otherwise 
uniform medium, thus dealing just with the situation typical for sonolumi- 
nescence.) Formally, in the presence of curved boundaries one is confronted 
with divergences when summing over all angular momenta up to infinity. This 
kind of divergence is usually absent when one deals with plane boundaries. 
Physically, the divergences are coming from the fact that phenomenological 
electrodynamics, implying use of the permittivity concept, becomes inappro- 
priate at small distances. There exists a natural cutoff in the material, of 
the order of the intermolecular spacing, and in practice some kind of reg- 
ularization has to be invoked in order to deal with the divergences in the 
formalism. More accurately this cutoff is the molecular diameter, as is most 
easily seen from the statistical mechanical method. By means of this method 
the electromagnetic Green function can be associated with the pair correla- 
tion function between dipole moments. The latter quantity is zero inside the 
hard cores, and its deviation from from the "ideal" Green function extends 
typically a few molecular diameters outwards from the molecule. 

In the case of nondispersive media, the use of zeta-function methods has 
proved to be very useful for force, or energy, calculations. The field theory 
approach to the Casimir problem has been considered at various places; in 
addition to the references above we may mention Refs. [4-16]. (This list is not 
intended to be complete; it covers mostly treatments of nonmagnetic media, 
and does not include the bulk of papers devoted to studies of the special case 
of media that satisfy the condition e\i = 1. A very extensive list of references 
is given in the recent report of Nesterenko et al. The Casimir energy 



E calculated by field theoretical methods at zero temperature for a dilute 
nondispersive compact sphere of radius a and refractive index n is positive, 

E= 23_hc ( 1)2 (1) 
384vra ; ' V ; 

corresponding to an outward force. 

Instead of making use of field theoretical methods for continuous matter, 
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one may alternatively use quantum statistical mechanical methods. We shall 
in the first sections below consider methods that were developed by H0ye and 
Stell, and others. Basic references to this kind of theory are ]T7J and [|Tg|]. 
In the Casimir context, H0ye and Brevik |i| used the quantum statistical 
mechanical path integral method to calculate the van der Waals force between 
dielectric plane plates. Recently, we have applied the same method to a 
single compact spherical ball pOfl . This statistical method, although probably 
not so well known as the field theoretical methods, turns out to be quite 
powerful. Thus, we can use it to calculate explicitly the short range terms in 
the single sphere's free energy, and verify how the repulsive Casimir surface 
force as calculated by field theoretical methods is simply a residual, cutoff 
independent, term in a complicated expression containing many terms. Cf. 
in this context also Refs. f2lfl , [|16|1 , and 0. 



It is now natural to ask: what is the experimental status in this field? 
Recently, there has been an impressive improvement of the experimental 
accuracy as regards force measurements; one has been able to verify the 
theoretically predicted Casimir forces, lying in the piconewton range, up to an 
accuracy of about 1 per cent. In Ref. |E| the Casimir force was demonstrated 



between metallic surfaces of a sphere above a disk using a torsion pendulum, 



whereas in Refs. Ifffil , [24| an atomic force microscope was used. 



One important lesson is to be learned from these experimental works is 
the following: they alway involve two (in principle there may even be more) 
bodies. The Casimir surface force on a single sphere is not measured. There 
seems not even to be an idea of how to measure such a force; probably this 
reflects simply that the force concept as such is not observationally well de- 
fined. Thus, in order to keep contact with experiments, at least in principle, 
one ought to consider at least two bodies. And this brings us to the theme of 
the present paper, namely to calculate the mutual free energy for a system of 
two spherically-shaped concentric nonmagnetic dielectrics. We will envisage 
that there is one compact sphere for r < a, and one semi-infinite similar 
medium for r > b, so that there is a vacuum gap of width d = b — a in be- 
tween. There will be an attractive Casimir force between the two media. One 
may object that there is still no straightforward way to imagine measuring 
such a force; however this does not create difficulties for our main purpose, 
which is to calculate the Casimir force in a setting which maintains spherical 
symmetry and yet avoids the complications with internal, cutoff dependent, 
forces. 
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In the following four sections we shall deal with the quantum statistical 
mechanical theory, with an emphasis on the static limit (zero Matsubara 
frequency, where derivations are more simple). Thereafter, for comparison we 
consider the alternative field theoretical approach, limiting us for simplicity 
to the case of perfect metallic walls at r = a, b. The resulting expressions 
for the free energy are Eq. (18) for the static case (Matsubara frequency 
equal to zero) and Eq. (40) for the time-dependent case (general Matsubara 
frequency). These expressions are obtained within the statistical mechanical 
approach. Within the field theoretical approach, the finite-temperature free 
energy is given by Eq. (68) assuming, as mentioned, perfectly conducting 
walls. 

There is one notable difference between the statistical mechanical ap- 
proach and the field theoretical approach, as far as the free energy is con- 
cerned. In the first of these cases the method is basically more simple, at 
least in principle, as one needs only knowledge about the mode eigenvalues of 
the oscillating dipole moments in the dielectric medium; cf. Eq. (3). These 
eigenvalues are again related to the pair correlation function. In the second 
case the calculation is more indirect, as one first calculates the surface force 
(arising from the mutual interaction) on the outer surface implying use of 
Maxwell's stress tensor, and thereafter relates the force to the free energy via 
integration of Eq. (67). That is, the field theoretical method involves use of 
the two-point functions for the electric and magnetic fields. As we also want 
to show that the results obtained by these widely different methods are in 
agreement, we consider some cases that are easy to analyse analytically, by 
both methods. 

We ought to stress again the conceptual difference between the two meth- 
ods studied in this paper. By the field theoretical method the Casimir effect 
is regarded as the energy shift due to the frequency eigenvalues of the quan- 
tized electromagnetic field in the presence of dielectric media. By the more 
recent statistical mechanical method this energy shift is regarded as a con- 
sequence of the dipole- dipole interaction between oscillating dipole moments 
embedded in the molecules of the media. The latter viewpoint can be real- 
ized by exploiting the analogue that exists between phonons in a solid, and 
electromagnetic waves or photons in vacuum. Then impurities in the solid 
will be the analogue to dielectric particles. These impurities will couple to 
the phonons of the solid and modify their frequencies. If, however, one wants 
to do the statistical mechanics of the latter system, e. g. to calculate its free 
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energy, then one can first integrate out the coordinates of the pure solid that 
appear in the path integral representation of the quantized problem. With 
harmonic oscillators one encounters in this way Gaussian integrals which can 
easily be calculated. The result is independent of the impurities, except that 
interactions between them are introduced. Thus, the resulting change in free 
energy can be related to these induced induced interactions (being in general 
time-dependent) in the system of impurities. Likewise, in this picture the 
dipole-dipole interactions are related to the Casimir free energy. 
We employ Gaussian electromagnetic units in this paper. 

2 General remarks 

Consider the free energy F(T) due to the mutual interaction between two 
spherical dielectric bodies with concentric surfaces at r = a and r = b. The 
attractive Casimir force between the surfaces, per unit area at the outer 
surface, is equal to / = —l/(A7fb 2 )dF/db. As shown earlier for the case of 
plates |19| , this Casimir force can be interpreted as the dispersion force aris- 
ing from thermal fluctuations of molecular dipole moments. By our quantum 
statistical mechanical considerations this also incorporates the quantum fluc- 
tuations at T = 0. That means, all fluctuations can be regarded as purely 
thermal for any T. The difference between classical and quantum situations 
is that in the former case these fluctuations vanish at T = while in the 
latter case they remain finite. In |^U| we considered the low density (or small 
e — 1) version of the single-body problem showing, as mentioned above, that 
the divergences are due to the continuum model of the medium. A cutoff 
in length scale is needed. For realistic systems this cutoff is determined by 
the molecular hard cores through their influence upon the pair correlation 
function. For two polarizable particles the free energy due to their mutual 
attraction can be written as |L9], 

1 00 1 1 
PF = --Y, -(«i^« 2 ^) n = - ln(l - a^arf), (2) 

Z n=l n 1 

with (3 = 1/ksT. The ip represents the potential energy of the dipole-dipole 
interaction which for general K (see Eq. (3) below) is given by Eqs. (6) and 
(7) below. Now the two particles can be generalized to and regarded to be 
our two spherical bodies, in the same way as two semi-infinite parallel plates 
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were treated in [ I9[ . That means, the two spherical bodies are regarded as 
two particles with many internal degrees of freedom. In this way Eq. (2) 
becomes a short hand notation wherein ip represents the interaction between 
two points in the two bodies (over which we integrate), and the polarizabili- 
ties cti and a 2 become the respective internal correlation functions of the two 
bodies with their mutual interaction ip switched off. As noted in |19[ the ex- 
pression (2) is formally exact for coupled harmonic oscillators, i.e., the model 
that we are employing for the polarizable particles. In terms of graphs, the 
expression (2) represents the ring graphs in the 7— ordering for the long-range 
forces, 7 being the inverse range of interaction [Q. For coupled oscillators 
Eq. (2) above and Eq. (3) below are exact results [j2~7j] . 

Extending to the quantum mechanical case, Eq. (2) generalizes to 

(3F = - ^2 ln (l - olxk^kolzk^k) , (3) 

1 K 

where K = 2im/ (3 with n integer (i.e., n e (—00, 00)). Note that K = h( n , 
where ( n is the Matsubara frequency ( n = —iuj and u is the frequency. 

For low density (or small a) only the first term in the sum (2) is needed. 

pOj and found there, after some trans- 
The first term means 



This is the situation considered in 
formations, to be in agreement with earlier works 
simply that one takes the (radiating) dipole interaction squared, average (in- 
tegrate) over the fluctuating dipole moments, and finally integrate over the 
two media. Thus 



F = p 2 dridr 2 



where p is the particle density, r% < a, r 2 > b, and [JIS 
0* = 42X Wl K (r)+^l K {r) 



A 



(4) 



(5) 



where r = r 2 — rj.. (Here ri and r 2 are positions in different media, so double 
counting does not occur.) The radiating dipole interactions used in (5) can 
be written as 



with 



V(12) = ip DK (r)D K (12) + il>AK(r)A K (12), 



(6) 



D K (12) 
A K (12) 



3(ra 1K )(r d 2K 



O'lK a 2K, 
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Here the hats denote unit vectors, and is the Fourier transform of the fluc- 
tuating dipole moment of particle number i in imaginary time; cf. Eq. (5.2) 
2§. Explicitly, from Eq. (5.10) in |25J, 



m 



^ DK {r) = - e —{l + T+ l -T 2 
Wr) = -^|r 2 + ^(r), (7) 



with 



iur Kr 



for 9f(w) < (or K < 0), 



c he 

and —K — > |if| when extending to K > in (7) (see Eq. (5.11) in f25[). 

For general we are not able to calculate the integral (4) in a direct way 
(but we can calculate it indirectly for arbitrary density, as will be argued 
later, in Sect. 4). We can integrate, however, for K = 0. The integral of 
interest then becomes 

1= [ dridr 2 —, r 2 = r\ + r\ — 2r\r 2 cos 9. (8) 

Jri<a,r2>b T 

This integral can be evaluated in closed form (second reference). How- 
ever, in Sect. 3 we want to extend this low density evaluation to the case of 
arbitrary density or arbitrary e. To do so, we need the contributions related 
to the various spherical harmonics. Thus we will here perform an expansion 
of the integral. This will also be used as an independent verification at the 
end of Sect. 4, where the more general theory developed there for arbitrary 
values of K turns out to yield the correct result when K — ► 0. 

Using spherical coordinates to integrate over the angle 9 between ri and 
r 2 we first obtain (x = — cos 9) 



J 



dx 1 



i (r\ + r\ + 2rir 2 x) 3 4rir 2 



1 



(r 2 - r x ) 4 (r 2 + r x ) 4 



1 g (2/ + 2)(2Z + l)2Z fr_i\ 2l -\ {9) 



performing a series expansion to make it easy to relate to the result for high 
density. Then, 

^2 



37T 2 



/ ridr x I ridr 2 J = > ^- -a h (10) 

Jo 1 Jb 2 3 ^ 21 + 1 " v ' 
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with a x = {a/b) 2l+1 . 

Inserted in (4) we obtain (a = oc) 



I3F = f3p 2 -a 2 -2-1 = ~(e - l) 2 V ^-^oj. (11) 
^ 2 2 v ; fr[ 21 + 1 v ; 

Note that for small e — 1, Eq. (11) will be the high temperature result for 
which only = contributes. This high temperature result at low density 
may in itself be of limited interest as it does not go beyond earlier results. 
But here we use it as a basis to make further developments. So in the next 
section the formalism is generalized to arbitrary density or e, although it is 
still restricted to K = 0. In Sect. 4, a derivation that encompasses both 
arbitrary e and K is given. 



3 The static case 

For simplicity we first consider the static case, by which we mean that the fre- 
quency is zero (K = 0). Then the electromagnetic dipole-dipole interaction 
is the well-known static, time-independent (also called instantaneous) one. 
By the time- dependent case we mean the general situation with K ^ 0. Then 
the dipole-dipole interaction will be the radiating or dynamical electromag- 
netic field where the time delay due to the finite speed of light is involved. 
Note that in general both the static and the dynamic cases contribute to 
the Casimir effect. The former, being proportional to T, contains the whole 
effect when T — > oo, but vanishes when T — ► 0. 

For general e one should sum up the series in Eq. (2). This will not be 
a simple task. However, one can include a strength factor A along with the 
perturbing interaction ip and differentiate (2) (or (3)) to obtain 

OF 

^qvU=i = -V>(aiV>cO!2), ( 12 ) 

where 

/ v> 

Vc 



1 — (XY^oii^ 

Here the axij) c a,2 will be the pair correlation function for the fluctuating 



dipole moments. As shown in Appendix A in |l9j the ip c (apart from a 
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simple factor) is the Green function for the electromagnetic problem with 
the dielectric medium present while ip is the one for vacuum. Thus we can 
utilize Maxwell's equations for electrostatics to obtain this zero frequency 
Green function or ip c in the presence of two dielectric spheres. 

The electrostatic potential $ fulfils the Laplace equation V 2 $ = with 
e =const. Splitting off the spherical harmonic factor Yi m = Yi m (9, <p), 



$ = $ l (r)Y lm (0,<p), 
we can write the radially dependent term in the form 



(13) 



r<a 

C^) l+1 + C,^)\ a<r<b 
[D(s) , b<r. 



(14) 



From the boundary conditions the coefficients can be determined. We give 
the coefficient D belonging to the exterior region 



D 



{21 + If 



At 



where 



( e -l)2(Z + l)Z (1-^(7,)' 
s-l) 2 (l + l)l 



and o\ 



2«+l 



(15) 



(16) 



[e(l + l) + l](el + l + l) 

The coefficient D represents the change of the field for r > b relative to 
the 5 = 1 case for a given point source. Via Eq. (12) the free energy is now 
obtained in a straightforward way. For small e — 1 the quantity (12) is twice 
the quantity (11). Thus the expression (12) becomes 



8 dF \ - Y7 2 f i I) Am 



(17) 



As we will argue below Ai<ii will be proportional to the strength factor 
squared, A 2 . Thus integrating (17) we obtain the free energy 



1 oo 

/3F = -£(2Z + l)ln(l-^), 

2 i=i 



(18) 
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which clearly yields an attractive force between the two spherical bodies. 

It should be noted that the series (18) is convergent as it contains no 
self-energy. From Eq. (16) it is seen that A\ < 1 and o\ decays like an 
exponential as I — > oo, implying AiO\ < 1. Further, no cutoff is needed as 
Eq. (18) only includes the free energy shift due to the mutual interaction 
between the two bodies. A "self-energy" would arise if interactions within 
each body were considered. But yet the necessary minimum distance between 
molecules would prevent the latter expression from diverging f20f. 



4 Further analysis of the static case 

The time-dependent case (K ^ 0) will be more complex to handle as we 
are not able to perform analytically the generalization of the integration (9) 
that gave the results (11) and (18). We find however that this case can be 
handled indirectly, noting that the quantity 

M = a x i\)a 2 i\) (19) 

in Eqs. (2) and (12) can be regarded as a matrix. We want the trace of these 
expressions (as well as the expression (3)), which amounts to integrating 
over positions and dipolar moments of the particles. The matrix can be 
transformed into a diagonal matrix A through some matrix S, 

M = SAS- 1 . (20) 

Then 

Tr(M q ) = Tr(A q ) = ^Af, 

i 

where Aj are the diagonal elements of A. Also, 

(M% = J2Su(S-%Xl (21) 
i 

Thus to obtain the free energy (2) or (3) we only need the eigenvalues A/. 
Use of the spherical harmonics Yi m for our present problem produces such 
a diagonalization and, as the results (17) and (18) show, the AiO\ represent 
these eigenvalues. The prefactor 21 + 1 is simply the degeneracy factor. 

However, without performing the integration (10) the identification of 
Ai<j\ with the appropriate eigenvalues is not obvious and cannot be concluded 
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from Eq. (15) alone. Then we turn to Eq. (12) and consider the correlation 
function (or the equivalent Green function) which can be expanded as 

«i^ c «2 = : aiij °! 2 . = a^a 2 £ M n (22) 

1 - (XxtyOLilp ^ 

with M = a\ipa 2 ip. Applying S the M is made diagonal such that 

S^ 1 aiip c a 2 S = S^ 1 aiipa 2 S- — — , (23) 
and A or its eigenvalues Aj can be identified via the ratio 

5 - 1 «W^2 s= 1 (24) 

a\ipa 2 1 — A 

Here the D in Eq. (15) represents the numerator (the full correlation func- 
tion), while the denominator can be identified with the first term in a chain 
bond expansion with one single potential bond ip and two hypervertices a% 
and a 2 (or correlation functions for the two media with their mutual in- 
teraction switched off). By chain bond expansion we mean the graphical 
representation of the terms in the expansion (22), where hypervertices aij al- 
ternate with potential bonds ijj. This notation has its origin in the statistical 



mechanical theory of fluids ||26|| . To go beyond standard mean field theory 
the Mayer graph expansion can be rearranged, and the chain bond will then 
become the leading correction (for forces of long range) to the correlation 
function of the reference system (e.g., hard spheres). In the present case 
with harmonic oscillators the expansion turns out to be exact for the pair 
correlations of amplitudes (corresponding to Gaussian fluctuations). 

The first term of the expansion (22) is now obtained by considering the 
two spherical bodies separately, or equivalently by considering Eq. (14) first 
with a = and thereafter with b = oo. Then there will be no multiple ijj 
bonds going back and forth, as there are no longer two media present. First 
take away the inner sphere, which is done by putting a = 0. Solving for D 
one obtains (B = 0) 

D — Dq — -Ji^C. (25) 

The amplitude ratio D /C will represent a 2 (or a\). Secondly, take away 
the outer sphere, which is done by putting b = oo. Solving for C one then 
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obtains (-0 = 0, d = 0) 



c = a 



21 + 1 



(26) 



oo 



el + 1 + 1 



The amplitude ratio 0^/(1/ e) (see Eq. (14)) will represent a\ (or a 2 )- 
Thus D Q will represent aiip(X2 as if) is represented by 1/e. Likewise the D as 
given by Eq. (15) represents the full correlation (or Green) function as both 
spheres are present, by which the ratio D/D yields the sought eigenvalues 
of Eq. (24). 

Combining Eqs. (25) and (26), 



Relating this to Eqs. (15) and (24) we see that the eigenvalues of A are 
\i = Ai<Ti, D/D = 1/(1 — A/). Thus we recover the results (17) and (18) 
when A; is used in the expressions (12) and (2), taking into account the 
degeneracy factor 21 + 1. That means, the present indirect approach is fully 
consistent with the explicit integration (10) that led to the same results. 

5 The time-dependent case 

Including the time dependence the solutions of Maxwell's equations become 
more complex. One has to solve the vector wave equation, and the fields have 
to satisfy the boundary conditions at the two surfaces. Again the spherical 
harmonics Y[ m can be used, and the remaining problem how to decompose the 
vector fields parallel and transverse to the spherical surfaces is conveniently 
dealt with in terms of the TM (transverse magnetic) and TE (transverse 
electric) mode Application of the angular momentum operator L = 

(l/«)r x V (with h = 1) creates a vector normal to r, i. e. r • L = 0, and is 
thus parallel to the spherical surfaces. As L commutes with the V 2 operator 
of the wave equation, and as L does not contain differentiation with respect 
to r, the wave equation has TM solutions of the form 



D, 



(21 + 1) 2 



(27) 



'o = 



[e(l + 1) + l](el + I + 1)' 



B 



$(r) 



(28) 



r 
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where $(r) is some function of r. Likewise the TE solutions follow with B 
replaced by E. 

The E field is now obtained from 

£ <9E 

V x B = -— = -ikeE, (29) 
c dt K J 

where k — u/c. Thus we need the formula p8| 

d 

iVxL = rV 2 — V(l + r— ). (30) 

or 

Now applying boundary conditions on the spherical surfaces, we find that 
the condition on the radial component of E coincides with that of B, so that 
we need its component Ej_ transverse to r which comes from the last term 
in (30) where only derivatives with respect to the polar angles are needed 
from the V operator. The latter again act only on Yi m which are the same on 
both sides of the interfaces and can thus be disregarded as far as boundary 
conditions are concerned. Therefore we are left with the r-dependence of E 
, where the term of interest is given by 

< 1+ 4) ( ^)=*m. (3i) 

The solutions of the wave equation for a given frequency are the Riccati- 
Bessel functions. As independent pair of functions it is convenient to choose 
the functions that are proportional to rji(kr) and to rh!^\kr)\ the first one 
because of its fmiteness at the origin, the second because of outgoing bound- 
ary conditions at infinity. After frequency rotation, and convenient nor- 
malization, these are the functions denoted by si and e\ in the field theory 
section below, in Eq. (53). For simplicity we will in the present section omit 
the subscript I. We will let subscripts a, b refer to functions taken at r = a, b, 
and add an extra subscript e to indicate that the function is taken inside a 
dielectric medium. Like Eq. (14) we can now write 




$(r) = { Ce + ds, a<r<b (32) 



As compared to Eq. (14) the coefficient 1/e for r < a has been dropped since 
the $(r) represents the magnetic field, but has no further consequence as it 
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only affects the other coefficients by a proportionality factor e. Requiring 
continuity of the tangential components Bj_ and Ej_ across the surfaces we 
obtain the equations 

e a£ + Bs a£ = Ce a + ds a , 
- £ (e' a£ + Bs'J = Ce' a + C lS ' a , 
Ce b + Cis b = De b£ , 

Ce' b + C lS ' b = je' b£ , (33) 

where we emphasize that the primes here mean differentiation with respect 
to r. 

To obtain the eigenvalues of interest we now proceed as in Sect. 4. So 
like Eq. (25) we find from the two last members of Eq. (33) (i. e., a = 0) 

D = D = c 2 C, with c 2 = s 6bS \' e '^ b , (34) 

ee be s b - e b£ s b 

and like Eq. (26) we find from the two first members of Eq. (33) (i. e., b — oo 
and Ci = 0) 

C = Cod = Cu with Ci = e' a£ s at -e a£ s^ (35) 

Combining these we obtain, like Eq. (27), 

Do = cic 2 . (36) 

As explained below Eq. (24) and in connection with Eqs. (25) and (26), this 
.Do represents the chain with a single potential bond ip. The full chain bond 
(see explanation below Eq.(24)) will be obtained by solving Eqs. (33) as they 
stand. This yields 

D = (37) 
where / f f / 

^ _ ( £s 'a s ae ~ s a s 'ae)( £e 'b e be ~ e b e 'be) /gg\ 

d ( £e 'a s as - e a s' a£ )(ee b£ s' b - e' b£ s b ) 

are the eigenvalues of interest in the construction of the free energy; cf. the 
argument above Eq. (27) in Sect. 4. 
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The static case (uj = 0) is recovered by putting e e = e = l/r l and 
s E = s = r l+1 , which yields X £ i = Ai<j\ in accordance with Eq. (15). Note 
that Eqs. (33) are somewhat different from those used in the static case as 
e is replaced by 1/e while / and I + 1 are interchanged, but the result is the 
same. 

When u ^ there is also another set of modes, namely the TE modes. 
They are obtained by replacing B with E in Eq. (28), and by interchanging 
B and E in Eq. (29) while removing the factor e and the minus sign on 
the right hand side. Again imposing boundary conditions, Eqs. (33) are 
recovered, except that the factor 1/e is no longer present. Solving for D 
we recover the results (34)- (38) also, except that all factors e are no longer 
present. The eigenvalues of interest now become 

^ _ ( S a S a.£ ~ s a s ae)( e b e be ~ e b^be) {39) 
( e 'a S ae ~ e a,s' a£ ) (eb £ s' b — e' be Sb) 

With the eigenvalues (38) and (39), the expression (18) for the free energy 
can be extended in a straightforward way to the time-dependent case, and 
we get 

1 oo 

P p = o E £(2' + 1)N1 - U + M 1 - Ai)], (40) 

Z K 1=1 

where the prefactor 21 + 1 is again the degeneracy factor, and K = 2im/f3 
with n integer (i.e., n e (— oo, oo)). 

Finally it can be noted that for two parallel plates separated by a distance 
d the free energy can be written in a similar general form. This energy can 
be found by integrating the surface force given by Eq. (2.9) of Ref. [fl"9~ |. This 



surface force is the famous Lifshitz result [2?J. The parallel plates result for 



the mutual free energy per unit area becomes 

K 



1 r°° 

PF=—Y: / Ml - K q ) + ln(l - X q )] q dq, (41) 

47T K JCn/c 



where 
with 
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1/2 

, K = —iTiuj = h( n = 2nn//3. 

It should be noted that Eqs. (40) and (41) are general results, valid for 
arbitrary permittivity e(u) and temperature, and they give the free energy 
due to the mutual interaction between the two media. There are no diverging 
self-energy terms. From Eq. (37) one must expect X e i, Xi < 1, and in the K = 
case (18) one has the converging factor o\ = (a/b) 2l+1 . This convergence 
will not be weakened when K ^ as then exponential factors also enter the 
solution of the wave equation for imaginary frequency. 




6 Field-theoretical approach: the surface force 

We now consider, as an alternative, the field theoretical approach to the 
same physical system, with the simplification, however, that the compact 
media are perfect conductors (e = oo). Thus we can compare with the above 
general result in this special case. We shall make use of the local Green- 
function method, as developed in particular by Schwinger and his school. A 
basic reference to this kind of theory applied to the case of spherical sym- 
metry (a perfecly conducting shell) is Milton et al. To our knowledge 
Milton was also the first to apply this theory to the compact ball problem 
J|]. Generalization of the theory, so as to take into account electrostriction, 



was made by Brevik [31]. Later references are |J-|| and [10, |TT[. (This list 
does not include the main part of the references dealing with efi = I media, 
as well as papers dealing with the mode summation method.) We now put 
h — c — ks — 1. 

Once the assumption about perfect conductors is accepted, the formalism 
becomes relatively simple. Since all fields in the regions r < a and r > b 
are equal to zero, we have to consider the fields in the vacuum gap only. 
The Green function F(x,x') for two spacetime points x and x' has a Fourier 
transform r(r, r', u) defined by 

r(x,*')= r ^V^r(r,r>), (42) 

J-oo 2,71 

with r = t — f. Note that the convention of Fourier transform used here 
implies a change of sign of to (i. e., u — > —uj), compared to the definition 
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used in the preceding sections, e. g. Eq. (7). The governing equation for F, 
as following from Maxwell's equations, is 

V x V x r(r, r', w) - co 2 T(r, r', u) = uj 2 15(r - r'), (43) 

and the spectral two-point function for the electric field components is 

i{E l {v)E k {v')) u) = (47r)r ife (r,r» (44) 

(the prefactor An appearing because of our present use of the Gaussian 
system of units). There are two scalar Green functions, F t (r, r') and G t (r, r'), 
since there are two independent field modes. The connection between these 
functions and the spectral two-point functions is 



t(E r (r)E r (r')) u 



(4 



21 + 1 



rr' " An 



(45) 



2Z + 1, 



i(E ± (r)E ± (r')) w = (4vr) £ ^-[^(r, r' 



i=i 



hir r t> r ' Gi{r ' r% (46) 



=1 



47T 



(47) 



i(H ± (r)H ± (r')) w = (Air) f> 2 Q(r, r') + ±-^ r -^r' F^r')}. (48) 



Here, it is assumed that the vectors r and r' lie in the same angular direction. 
The radial difference r — r', however, does not have to be small. 

For simplicity we shall denote the scalar Green functions generically by 
Aj(r, r') (thus A; is either F t or Gi). Their governing equation is 



d 2 2 d 2 1(1 + 1) 



A,(r,r') = -^5{r-r'). (49) 



The solution contains spherical Bessel functions, and has the general form 



A,(r,r') 



ik 

1-CtC 



I^II 



j,(Ar<) - C / / i [ 1) (A;r < ) ^(fcr>) - C^(^>) 



(50) 
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where k = Ci and Cji being constants. This form satisfies the discon- 
tinuity condition following from Eq. (49) on the radial derivative of A/ at 
r = r', with the Wronskian W{ji(x),h\ 1 \x)} = i/x 2 . Taking into account 
the boundary conditions at r = a, b we can determine the constants: for the 
F (or TE) mode we get 

~ sijka) =, . . ei(kb) 

CiF{ka) = Kka)' ClIF{kb) = iwr (5 ) 

whereas for the G (or TM) mode 

=S- e -^=m- (52) 

Here, we let prime mean differentiation with respect to the whole argument. 
We now perform a complex frequency rotation, uo — > iu, k — > i\ou\ = ik, 
and replace the conventional Riccati-Bessel functions §i(x) = xji(x), ei(x) = 
xh\ 1 \x) by new ones s;,e/ defined according to 



si(x) = (-i) l+1 si(ix) = ^I u (x), 




e l (x)=i l+1 e l (ix) = J — K u (x). (53) 



Here v — 1+1/2, I v and K v are modified Bessel functions, and the Wronskian 
of importance now is W{si ) e{\ = —1. The frequency rotation implies that 
we replace the "tilde" constants C by new constants C, in accordance with 
the relation C(ix) = (—l) l+1 C(x). Explicitly, 

c '^ = Wr c " a(v) = WY (54) 

where we here and henceforth let x and y be defined by x = ka, y = kb. 

We now return to the expansions (45)- (48) for the spectral two-point 
functions. Of main interest for us are these functions when the points r 
and r' are close to each other, but not overlapping. We moreover set the 
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time-splitting parameter r — t — t' equal to zero. Substituting the two- 
point functions in Maxwell's stress tensor we can calculate the surface force 
density on either of the two surfaces. We choose the outer surface r = b, 
since it will then become easy to relate the force to the free energy. Writing 
for simplicity (E^(r)) instead of (E r (r)E r (r')) r r^ r , we obtain for T = the 
various two-point functions at r = b— : 

lE 2 (h u _ ( 47 r dy^2l + 1 Sl (y) - C IG (x) ei (y) 

(El(b-)) = (H?(b-)) = 0, (56) 

(H 2 (b-\\ ~t 4?r ) d V 21 + 1 \ Sl ^ - ClG ^ 61 ^ i " C IF (x)e\(y) 
{ ±{ ]1 ~ ^ Jo Vy ti 4vr ^s'^-dG^y) ^s l {y)-C IF {x)e l { V y 

(57) 

(the prefactors 4n again reflecting the Gaussian units). 

Using Maxwell's stress tensor we can write the surface force density on 
the outer surface as 

f b = -^(E*(b-)) + ±(Hi(b-)). (58) 

Substituting Eqs. (55) and (57) into Eq. (58) we obtain, when taking into 
account the governing equation for the Riccati-Bessel functions, 

s' l '(y) = (l+l(l + l)/y')s l (y) (59) 

(and similarly for ei(y)), that 

-1 r ^2l + l l s' l (y)-C IF (x)e' l (y) , s'{(y) - C IG (x)e'{(y) } 



f - 1 r v , v ^ ^ + y s Ay)-^iF(x)e' l (y) sj[ 
Jb 2nb*Jo gy t[ 4vr ^^-CMxMy) 



4tt s t (y) - C IF (x)ei(y) s' t (y) - C /G (x)e',(y) 

(60) 

From this expression it is apparent how both modes F and G contribute to 
the force. 

We can write Eq. (60) as a sum of the following two terms: 

h = if + ft\ (6i) 
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where 



^_^L r d k^ 2 l±l! i nm _^)^M in^HMn f63l 
h -2nVJo 4n db ln{[1 ei (x) Sl (y) l[1 e[{x) aj(y) J *' (63) 

with y = kb (the operator d/db is taken at constant value of a). The 
expression (62) is the same as the inner contribution to the surface force on 
a perfectly conducting shell |30], [5|. This term does not involve the interaction 



between the two media, and will be discarded in the following. Of interest 
for us is the interaction term (63). As si(y) = \e v and ei(y) = e~ y for large 
y, it is evident from (63) that /& — > if the outer surface recedes to infinity 
while the inner surface is kept constant. This is physically as it should be. 

A remark is here in order, concerning the physical meaning of the two 
force terms in Eq. (61), in particular the possibility of making measurements. 
It ought to be emphasized, first of all, that the term fjf^ is a mathemati- 
cal construct. It does not seem to be possible to measure this term, not 
even in principle. If one imagines the case of a perfectly conducting sin- 
gular shell with radius a, the case studied in f3(| and also in ||, then 



has to be taken together with a similar term on the outside to make up a 
surface force that in principle ought to be measurable. However, as far as 
we know, no measurement has so far been made for even such a complete 
shell. In our case, what is measurable, at least in principle, is the interaction 
term (63). For instance, one might envisage to measure the attractive force 
between a micrometer-sized conducting sphere and a semi-spherical trough 
in a conducting plate, thus some kind of generalization of the atomic force 



microscope measurements reported in [23] and [24 



The expressions above refer to zero temperature. The transition to finite 
temperatures is made by means of a discretization of the frequencies, 

k -> K = 2<im/(3, x -> Ka, (64) 

with n an integer. The rule for going from frequency integral to a sum over 
Matsubara frequencies is 

P n=0 



20 



where the prime on the summation sign means that the n = term is 
taken with half weight. The finite-temperature force expression accordingly 
becomes (recall that v — I + 1/2) 



tint 

Jb — 



2nb 2 p 



-1 



OO OO Q 
71=0 1=1 UU 



si(x) ei(y) 
ei(x) si(y) 



][1 



s'i( x ) e 'i(y) 



]}, 



(66) 



where now x = 2-Kna/P, y = 2imb/ (3. 

As for measurements of a force like fl nt there are subtle problems although 
the free energy (Eq. (68) below) is well defined. One may imagine that the 
two spherical media are liquids, and that the outer spherical shell can move. 
(The inner shell can also move if liquid is added or removed through a small 
pipe.) However, there is an extra complication compared to the case of 
parallel plates, namely the change of radius of a spherical surface. This will 
change the free energy associated with the surface tension. Although the 
latter change of energy is finite, a more precise evaluation of it will obviously 
be a complex task. The molecular structure has to be taken into account at 
the surface where the density changes abruptly. 

It can be noted, however, that the sum of the two surface energies must be 
the opposite of the general result (40) for a — b, i. e., with the two spherical 
surfaces fused together. One is then left with no surfaces at all, and the 
bulk free energy applies everywhere. However, simply putting a = b in a 
continuum approach does not work; Eq. (40) will diverge unless a cutoff is 
introduced for large / to mimic the molecular diameter. Whether the self- 
force in Eq. (62) is possibly related to such a surface tension energy or 
not is a problem that we have so far not been able to clear up. 



7 The free energy, analytically and numeri- 
cally 



Since we have calculated the force density on the outer surface due to the 
mutual interaction, it is easy to derive the corresponding expression for the 
interaction free energy F mt . We imagine the outer surface to be displaced 
by a small amount db, while the inner surface is kept constant. The relation 



Jb 



■int 



I Q pint 



(67) 



4vr6 2 db 
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is integrated from b to infinity, noting that F mt = at b = oo. Making use 
of (66) we thus obtain, dropping the superscript 'int', 

valid at arbitrary temperatures. Comparing with the more general result 
of Eqs. (38)-(40) in Sect. 5 one finds that Eq. (68) agrees with these when 
e — > oo as then es' a s ae 3> s a s' ae , s' a s ae <C s a s' ae , etc. (with s a = si(x) etc. as 
introduced above Eq. (32)). 

At T = 0, where the free energy F is the same as the energy E, we can 
write 

B = ±r dx £, ln{[1 _£iM£iM ][1 _4M4M ]} . (69) 

iraJo JrJ e,(x)s,(y)" e',(x) s'fy)" ( ' 

Expressions (68) and (69) hold for arbitrary widths of the annular region. 
It is of interest, before embarking on numerical evaluations, to analyse some 
limiting cases by analytical means. The limiting case of immediate interest 
is that of a narrow slit, i.e., 

t= b — « 1. (70) 
a 

This case is motivated physically from the fact that the Casimir measure- 
ments are made for small separations only, and also because we have in this 
way the possibility to check our results against the standard results for par- 
allel plates in the limit when £ — > 0. 

At T = we find, when x and y lie close to each other || , 

si(x) ei(y) = s'^e'^y) = ^ 
ei(x) si(y) e[(x) s\{y) 

to 0(1/ v) in the uniform asymptotic (or Debye) expansion. Here 

21 + z 2 " Z ~u 



2£v / T + ^[l - + ...], z = - (72) 



Keeping only the first term, we find the T = interaction energy to be 

2 00 



E = —J2" 2 I d2ln(l-e- 2 ^). (73) 
na l=1 
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This expression can be processed further, if we make a power expansion of 
the logarithm and take into account the property Ya=i v 2 e~ v< ^ — > 2/0 3 when 
6^ [021. Then, 



E = -^-C(4) r d * = (74) 
2nae } Jo (1 + ^ 2 ) 3 / 2 180a e 3 ' 1 J 

which corresponds to the following interaction energy per unit surface (total 
surface area A = 4-na 2 , and d = b — a): 

E _ 7T 2 1 

A ~ "720 TP' (75) 

This simple calculation thus provides us with a satisfactory check: Eq. (75) 
is the conventional expression for the Casimir energy of two parallel plates. 
(The expression includes the effect of retardation. That is, the distance d is 
much larger than the characteristic wavelength of the absorption spectrum 
of the medium.) Thus, by retaining the first order term in £ we see that 
our theory reduces to the standard T = theory. Corrections to the theory 
arising from the curvature of the surfaces can in principle be worked out by 
going to larger powers in £. 

At finite temperatures, we obtain for the free energy 



(3F = 4£'X>ln(l - e ~ 2 ^ 2+n2 * 2 ), t = 2na/(3. (76) 

n=0 1=1 

This expression, as before, implies keeping of only the first term in the 
expansion (72), but it puts no restriction on the temperature. 

Let us consider the limiting case of high temperatures, first going back 
to the expression (68), holding for arbitrary widths d. For the highest tem- 
peratures (classical limit), only the lowest Matsubara frequency (n = 0) 
contributes. As x = nt, y = bnt/a, it is seen that we then need to evaluate 
Si and ei when the arguments tend to zero. As 



for small arguments, we have 

si(x)e t (y) = s'^x) e'^y) = ,a 2l+1 
ei(x) Si(y) e[(x) s[(y) b 



(78) 
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so that the contribution from n = becomes 



0F(n = 0) = £(2/ + 1) ln[l - Q) 2l+1 \. (79) 
z=i 

This is in agreement with our previus expression (18) (Ai = 1 when e — > oo), 
except from a factor 2. The physical reason for this is that both F and G 
modes contribute to Eq. (79), whereas only one mode contributes in (18). 
This artifact in Eq. (79) is related to the fact that e = oo is considered, while 
in Eq. (18) oj = is considered before the limit e — > oo is taken. 

In the case of a narrow slit we obtain from (76) the n = contribution 

oo 

PF(n = 0) = 2^>m(l - e~ 2 ^). (80) 
i=i 

Again making a power expansion of the logarithm, and observing the rela- 
tionship J2i=i ve~ l " t> = 1/0 2 when — > 0, we get 

m«=o) = -f^, (si) 



which corresponds to 

PF(n = 0) C(3) 1 



8tt rf 2 ' 



(82) 



Again, this is a satisfactory check, as Eq. (82) is the conventional high- 
temperature result for parallel plates. 

For a narrow slit we may also obtain the known result for a parallel plates 
configuration more generally. In the wave equation the term 1(1 + l)/r 2 ~ 
/(/ + l)/a 2 is replaceable with k\ where kj_ is the transverse wave vector, i.e., 

v 2 = (I + i) 2 ~ 1(1 + 1) = A; 2 a 2 . (83) 

When I is large we can regard it as continuous quantity, whereby the sum 
can be replaced by an integral. We have then 

YP l + x ) / ( 2/ + l ) dl = 2fl2 / k ± dk ±- ( 84 ) 

i=i J J 



24 



Further, 

27m , r ^ , 

P 



a 

or 



£v = —k±a = k±d, (85) 



2^v 2 + n 2 t 2 = 2qd, with q 2 = k\ + K 2 . (86) 
Insertion of this into Eq. (76) yields 

00 />oo 

f] F = 4a 2 J2' qdq ln(l - e~ 2qd ), (87) 

r> J Cn 



n=0 



using qdq = k±dk±. The surface force density thus becomes {b — > a) 

1 dF 2 00 r 00 e -2grf 
/b = "4rf¥ = "^S'l ^(l-e-^)' (88) 

with ( n = K = 2im/i3 being the Matsubara frequency. This is in agreement 
with Eq. (2.9) in |Tj| (A n = 5 n = 1 for e = 00). It is also in agreement with 
Eq. (3.8) in |33| (it should be mentioned that q 2 in our present notation, 
Eq. (86), is the same as k 2 in [[33], and also that the distance d above is the 
same as a in ]33| and (T^j). Note that whereas the expression (88) presupposes 
a narrow slit (large /), there is no restriction on the temperature. 

After these preliminary analytic considerations we now present numer- 
ically calculated results for the free energy, when the walls are conducting 
(e = 00). All results are given in nondimensional form. The numerics turned 
out to be rather demanding; as preliminary tests indicated that Matlab would 
be insufficient for our purpose we turned to standard FORTRAN routines 
and made use of them throughout. On a logarithmic plot with base 10, 
Fig.l shows how lg(— (3 Ft) = lg(— 2waF) varies with relative width d/a for 
various values of the nondimensional temperature t = 2ira/i3. At zero tem- 
perature, an integration routine was used for the x integral in Eq. (69). At 
finite temperatures, the double sum in Eq. (68) was calculated as it stands 
(thus without expansion procedures for the Riccati-Bessel functions), with 
use of the FORTRAN library for the Bessel functions of half-integer order. 
Allowing the numerical tolerances in the / sum as well as in the n sum to be 
equal to 10~ 6 , we found for the case of t = 0.01 and d/a = 0.1 the necessary 
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number of terms in the n sum to be about 415 000. For larger widths the 
necessary number of terms turned out to be considerably less; for instance 
for d/a — 1 and the same t the number was about 8900. 

A characteristic property of the curves in Fig.l is that they tend to overlap 
in the case of low temperatures. Thus the curve calculated for t — (via an 
x integral and an / sum) is indistinguishable from other curves calculated in 
the whole temperature region (via a double sum) up to t = 1. Numerically, 
for d/a — 0.1 the difference in lg(— (3 Ft) between the cases t — and t — 1 is 
found to be about 10~ 4 . For zero temperature, as mentioned, F is the same 
as the energy E. 

Whereas the representation in Fig.l is most useful for low temperatures, 
we show in Fig. 2 the representation of lg(— f3F) versus d/a. This is con- 
venient for high temperatures, since now the curves for high t stay inside 
the figure and tend to overlap. The curve calculated for t = 50 is actu- 
ally indistinguishable from the shown curve referring to t — 200. This fact 
reflects that F is proportional to t, thus in accordance with the behaviour 
of high temperature mutual free energies for classical harmonic oscillators 
in general. For a narrow slit, we expect that the calculation agrees with 
the approximate formula (82). We may check this in the case of t — 200, 
d/a = 0.05: the machine calculation then yields f3F = —249.7, whereas 
Eq. (82) yields (3F = -±C(3)(a/d) 2 = -240.4, thus an error of 4 %. The 
reason why the accuray here is only moderate, is most likely that the width 
parameter d/a = 0.05 is not small enough to represent a narrow slit to a high 
precision. Generally, it turned out to be difficult to calculate cases of higher 
temperatures or more narrow slits than those shown in the figures, without 
entering into special alterations of the FORTRAN routines. 

Finally, it is of interest to display explicitly the free energy's low-temperature 
plateau, and its high-temperature proportionality to t, for a fixed value of the 
relative width. This is done in Fig. 3, for the cases of d/a = {0.05, 0.075, 0.1}. 
The horisontal plateau is seen to prevail quite accurately until a gradual 
increase takes place in the region roughly between t — 10 and t = 30 
(lg30 = 1.477). It is notworthy that this behaviour is in agreement with 
the following simple physical argument: The most significant frequencies uj 
contributing to F are generally of the same order as the inverse width of the 
gap, i.e. uj ~ 1/d. Now, from Wien's displacement law we know that the 
maximum of a blackbody distribution occurs at a frequency of oo rn = 2.8//?. 
Temperature effects are expected to become significant when uj is comparable 
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to u rn . Putting uj ~ u m we obtain t = 2na/f3 ~ 2a/ d. In the present case 
this amounts to t ~ 20 — 40, (lgt ~ 1.3 — 1.6), which is seen to be in good 
agreement with the location of the shoulder in Fig. 3. Moreover, for higher 
temperatures, the proportionality of F to t is clear from the figure. 
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Figure Captions 



Figure 1 

Logarithm (base 10) of mutual nondimensional free energy, lg(— (3 Ft) = 
lg(— 2ttclF), versus relative width d/a for various values of the nondimensional 
temperature t = 2na/f3. For < t < 1 the curves are overlapping (i.e., only 
one curve is drawn), consistent with Fig.3. 

Figure 2 

Same as Fig.l, but lg(— j3F) is shown, as appropriate for the case of high 
temperatures. For t larger than about 50 the curves are overlapping (only 
one curve is drawn), consistent with Fig.3. 

Figure 3 

Variation of F versus t for the case of d/a = {0.05, 0.075, 0.1}. 
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Figure 1 , J.S. Hoye et al. 
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Figure 2, J.S. Hoye et al. 
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Figure 3, J.S.Hoye et al. 



